
;    s.ylib : routines for math functions
;
;    Copyright (C) 1993  Claus Vohwinkel
;
;    This program is free software; you can redistribute it and/or modify
;    it under the terms of the GNU General Public License as published by
;    the Free Software Foundation; either version 2 of the License , or
;    (at your option) any later version.
;
;    This program is distributed in the hope that it will be useful,
;    but WITHOUT ANY WARRANTY; without even the implied warranty of
;    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;    GNU General Public License for more details.
;
;    You should have received a copy of the GNU General Public License
;    along with this program; see the file COPYING.  If not, write to
;    the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
;
;    You may contact the author by:
;       e-mail:  vohwinkel@vxdesy.desy.de
;      us-mail:  Claus Vohwinkel
;                SCRI/FSU
;                400 Science Library
;                Tallahassee, FL 32306
;      


R0 RN 0
R1 RN 1
R2 RN 2
R3 RN 3
IP RN 12
SP RN 13
R14 RN 14
PC RN 15
F0 FN 0
F1 FN 1
F2 FN 2
F3 FN 3

        AREA |C$$code|, CODE, READONLY
 EXPORT yacos
 EXPORT yasin
 EXPORT yatan2
 EXPORT yatan
 EXPORT ysqrt
; EXPORT ycosh
; EXPORT ysinh
; EXPORT ytanh
 EXPORT ypow
 EXPORT yexp
 EXPORT ylog10
 EXPORT ylog
 EXPORT ycos
 EXPORT ysin
 EXPORT ytan
 EXPORT ymodf
 EXPORT yfmod
 EXPORT yceil
 EXPORT yfloor


 EXPORT exm1
 EXPORT twotox


; EXPORT yatanf
; EXPORT ysqrtf
; EXPORT twotoxf
; EXPORT yexpf


yceil STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 URDEP F0,F0
 NRME F0,F0
 mov pc, r14

yfloor STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 URDEM F0,F0
 NRME F0,F0
 mov pc, r14


yfmod STMFD SP!,{R0-R3}
 LDFD F1,[SP,#8]
 TEQ R2,#0
 MNFMIE F1,F1
 LDFD F0,[SP],#16
 DVFEZ F2,F0,F1
 URDEZ F2,F2
 NRME F2,F2
 MUFE F1,F1,F2
 SUFE F0,F0,F1
 mov pc, r14

ymodf STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 URDEZ F1,F0
 NRME F1,F1
 SUFE F0,F0,F1
 STFD F1,[R2]
 mov pc, r14


yacos STMFD SP!,{R14}
 BL yasin
 LDFE F1,PIO2
 SUFE F0,F1,F0
 LDMFD SP!,{PC}

ONE_O_SQT2S DCFS 0.707106781

yasin STMFD SP!,{R0,R1,R14}
 LDFD F0,[SP]
 MUFE F0,F0,F0
 RSFE F0,F0,#1
 BL ysqrtf
 LDFD F1,[SP]
 LDFS F2,ONE_O_SQT2S
 ABSE F3,F1
 CMF F3,F2
 BGT yasin1
 DVFE F0,F1,F0
 ADD SP,SP,#8
 LDMFD SP!,{R14}
 B yatanf
yasin1 DVFE F0,F0,F3
 BL yatanf
 LDMFD SP!,{R0,R1,R14}
 LDFE F1,PIO2
 SUFE F0,F1,F0
 TEQ R0,#0
 MNFMIE F0,F0
 mov pc, r14



 

yatan2 STMFD SP!,{R0-R3}
 LDFD F0,[SP],#8
 LDFD F1,[SP],#8
 BIC R0,R0,#&40000000
 ORR R0,R0,R2,LSR #1
 BIC R0,R0,#&20000000
 ABSE F2,F0
 ABSE F3,F1
 CMF F2,F3
 DVFLEE F0,F2,F3
 DVFGTE F0,F3,F2   ; REVERSE ORDER
 ORRGT R0,R0,#&20000000
 STMFD SP!,{R0,R14}
 BL yatanf
 LDMFD SP!,{R0,R14}
 TST R0,#&20000000
 LDFNEE F1,PIO2
 SUFNEE F0,F1,F0
 MOVS R0,R0,LSL #1
 LDFMIE F1,PI
 SUFMIE F0,F1,F0
 MNFCSE F0,F0
 mov pc, r14
 

PI & &00004000
   & &C90FDAA2
   & &2168C235


yatanf
 SFM F0,1,[SP,#-12]!
 LDMFD SP!,{R0-R2}
 TEQ R0,#0
 MNFMIE F0,F0
 BIC R2,R0,#&FF000000
 SUB R2,R2,#&3F00
 SUB R2,R2,#&FF
 MOVS R1,R1,LSL #2 ; WE SKIP THE FIRST DIGIT WHICH IS IMPLIED TO BE ONE.
 ADC R2,R2,R2
 B yatanx

; ------------------------

yatan
 STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8

 MOV R2,R0,LSL #1
 TEQ R0,#0
 MNFMIE F0,F0          ; WE MUST KEEP R0 FOR SIGN
 MOV R2,R2,LSR #20
 SUB R2,R2,#&FE
 SUB R2,R2,#&700
yatanx CMP R2,#7
 MOVGE R3,#24
 BGE BOK       ; WE ARE IN HIGHEST SECTOR
 RSB R2,R2,#6
 MOV R2,R2,LSL #2 ; * 4
 CMP R2,#32
 LDRLT R3,MASK1
 SUBGE R2,R2,#32
 LDRGE R3,MASK2
 MOV R3,R3,LSR R2
 TST R3,#8
 AND R3,R3,#7
 ADD R3,R3,R3,LSL #1 ; * 3
 BEQ BOK      ; THE BOUND IS OK
 ADR R2,COFFS1
 ADD R2,R2,R3,LSL #2
 lfm F1,1,[R2],#0
 CMF F0,F1
 ADDGT R3,R3,#3
BOK SUBS R3,R3,#3   ; WE MUST KEEP R3 FOR LATER ADDITION OF ATAN(X_I)
  BMI NOCORR
  ; R3 IS LOWER BOUND*3 = 0.. 3*7
     ; TO R3 BELONGS TAN(2*(R3/3+1)*PI/32)
 ADR R2,COFFS2
 ADD R2,R2,R3,LSL #2
 lfm F1,1,[R2],#0
 MUFE F2,F1,F1
 ADFE F0,F1,F0
 ADFE F2,F2,#1
 DVFE F0,F2,F0
 SUFE F0,F1,F0
NOCORR 
; WE HAVE THE NEW ARG IN F0
 ADR R1,ATCO6
 MUFE F1,F0,F0
 lfm F2,2,[R1],#24
 MOV R2,#4
 MUFE F2,F2,F1
 ADFE F2,F2,F3
ATA2
 MUFE F2,F2,F1
 lfm F3,1,[R1],#12
 SUBS R2,R2,#1
 ADFE F2,F2,F3
 BPL ATA2
 TEQ R3,#0
 ADRPL R2,COFFS3
 MUFE F0,F0,F2
 ADDPL R2,R2,R3,LSL #2
 LFMPL F1,1,[R2],#0
 ADFPLE F0,F0,F1
 TEQ R0,#0
 MNFMIE F0,F0
 mov pc, r14


ATCO6  & &00003FFB
  & &98FCA17B
  & &0B4B4F14
ZCO5  & &80003FFB
  & &BA1F194D
  & &4CBA52FF
ZCO4  & &00003FFB
  & &E38E1ED0
  & &E5D5B014
ZCO3  & &80003FFC
  & &92492486
  & &F9EEBEA5
ZCO2  & &00003FFC
  & &CCCCCCCC
  & &C8200EE7
ZCO1  & &80003FFD
  & &AAAAAAAA
  & &AAAA4DE2
ZCO0  & &00003FFE
  & &FFFFFFFF
  & &FFFFFFFF


MASK1 & &BCD6E77F    ; CONTAINS MASKS FOR EXP   -5  0   5  1  15  2  25  3
                     ;                           B  C   D  6   E  7   7  F
MASK2 & &0081192A    ;                         -45 -4 -35 -3 -25 -2 -15 -1
                     ;                           0  0   8  1   1  9   2  A
COFFS1             ;TAN( (2*I+1)*PI/32) I = 0..7
 & &00003FFB
  & &C9B5DC62
  & &D96D1310
 & &00003FFD
  & &9B5041AA
  & &E31EDE29
 & &00003FFE
  & &88D5B8C8
  & &41A75C35
 & &00003FFE
  & &D2180157
  & &21427E59
 & &00003FFF
  & &9BF7EC6C
  & &3B03EF3C
 & &00003FFF
  & &EF789DB9
  & &E0A03110
 & &00004000
  & &D2FACF48
  & &3006E648
 & &00004002
  & &A27362CA
  & &D7EE626A
COFFS2         ; 1/TAN(2*(I+1)*PI/32)  I = 0..7
 & &00004001
   & &A0DFF712
   & &123B8845
 & &00004000
   & &9A827999
   & &FCEF3242
 & &00003FFF
   & &BF90C712
   & &D3A30A83
 & &00003FFF
   & &80000000
   & &00000000
 & &00003FFE
   & &AB0DC155
   & &BFCC82F5
 & &00003FFD
   & &D413CCCF
   & &E7799212
 & &00003FFC
   & &CBAFAF02
   & &A98AC03E
 & &00000000
   & &00000000
   & &00000000
COFFS3         ; (I+1)*PI/16  I = 0..7
 & &00003FFC
   & &C90FDAA2
   & &2168C235
 & &00003FFD
   & &C90FDAA2
   & &2168C235
 & &00003FFE
   & &96CBE3F9
   & &990E91A8
 & &00003FFE
   & &C90FDAA2
   & &2168C235
 & &00003FFE
   & &FB53D14A
   & &A9C2F2C2
 & &00003FFF
   & &96CBE3F9
   & &990E91A8
 & &00003FFF
   & &AFEDDF4D
   & &DD3BA9EE
PIO2 & &00003FFF
   & &C90FDAA2
   & &2168C235


FSQINV TEQ R0,#&80000000
 TEQEQ R2,#0
 TEQEQ R3,#0
 ADD SP,SP,#12
 MOVEQS PC,R14
 TST R0,#&80000000
 MOVEQS R1,R0,LSL #2
 SQTNED F0,F0
 MOVNES PC,R14
 BIC R0,R0,#&C0000000
FD1 MOVS R3,R3,LSL #1
 ADCS R2,R2,R2
 SUB R0,R0,#1
 BPL FD1
 MOV R1,#&3F00
 ADD R1,R1,#&FF
 STMFD SP!,{R1,R2,R3}
 lfm F0,1,[SP],#12
 SUB R2,R0,R1
 B SQA

FZER1 TEQ R0,#0
 TEQEQ R2,#0
 TEQEQ R3,#0
 ADD SP,SP,#12
 MOVEQS PC,R14
 MNFE F0,#1      ; if we arrive here, something is screwed
 SQTD F0,F0
 mov pc, r14

ysqrtf sfm F0,1,[SP,#-12]!
 LDMIA SP,{R0,R2,R3}
 TST R0,#&C0000000
 BNE  FSQINV
 MOV R1,#&3F00
 ADD R1,R1,#&FF
 STMIA SP,{R1}
 TEQ R2,#0
 BPL FZER1
 lfm F0,1,[SP],#12
 SUB R2,R0,R1
 B SQA

ysqrt
 MOVS R2,R0,ASR #20
 BMI SQINVALID
 MOV R3,#&FF
 ORR R3,R3,#&300
 BEQ SQZER_OR_DEN
 BIC R0,R0,R2,LSL #20
 ORR R0,R0,R3,LSL #20
 EOR IP,R2,R3
 TEQ IP,#&400
 BEQ SQOVERFL_OR_NAN     ; THIS MEANS WE HAVE INFINITY OR A NAN ON INPUT
SQ1 STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 SUB R2,R2,R3      ; NEW EXPONENT CORRECTION * 2

SQA
; MOV R0,#2
 ADR R1,SQCO3
 LDFS F1,[R1],#4
;SQ2 ; unrolled - Mech
 FMLS F1,F1,F0
 LDFS F2,[R1],#4
; SUBS R0,R0,#1
 ADFE F1,F1,F2
; BPL SQ2
 FMLS F1,F1,F0
 LDFS F2,[R1],#4
; SUBS R0,R0,#1
 ADFE F1,F1,F2
; BPL SQ2
 FMLS F1,F1,F0
 LDFS F2,[R1],#4
; SUBS R0,R0,#1
 ADFE F1,F1,F2
; BPL SQ2

 MUFE F2,F1,F1 ; Y0^2
 MOVS R0,R2,ROR #1
 ADDMI R1,R1,#12
 ADFE F3,F2,F0 ; Z = Y0^2 + X
 MUFE F0,F0,F2 ; X*Y0^2
 LDMNEIA R1,{R0,R1,R3}
 MUFE F2,F3,#0.5 ; 0.5*Z
 ADDNE R0,R0,R2,ASR #1
 MUFE F1,F1,F3 ; Y0*Z
 STMNEFD SP!,{R0,R1,R3}
 MUFE F2,F2,F2 ; 1/4*Z^2
 LDFNEE F3,[SP],#12
 ADFE F0,F0,F2
 MUFNEE F0,F0,F3
 DVFE F0,F0,F1
 mov pc, r14


SQCO3 DCFS 0.0249374292634
SQCO2 DCFS -0.182747647015
SQCO1 DCFS 0.787923250524
SQCO0 DCFS 0.369956381214
SO & &00003FFF
   & &80000000
   & &00000000
SX & &00003FFF
   & &B504F333
   & &F9DE6484


SQZER_OR_DEN ;
 TEQ R0,#0
 TEQEQ R1,#0
 MVFEQSZ F0,#0
 MOVEQS PC,R14 ; input was zero
 MOV R2,#1
SQNORM MOVS R1,R1,LSL #1
 ADC R0,R0,R0
 SUB R2,R2,#1
 TST R0,#&100000
 BEQ SQNORM
SQNORM1
 ORR R0,R0,R3,LSL #20
 B SQ1

SQINVALID TEQ R0,#&80000000
  TEQEQ R1,#0
  MNFEQSZ F0,#0
  MOVEQS PC,R14
  B CALL_SUP

SQOVERFL_OR_NAN ORR R0,R0,#&40000000
CALL_SUP STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 SQTD F0,F0
 mov pc, r14


EXM1COFFS  ; P2 P1 P0 Q1 Q0
   & &00003FF9 ; 0.023094321272953855
   & &BD304D52
   & &E6A303D0
   & &00004003 ; 20.201700006953125
   & &A19D14E4
   & &ABBBB271
   & &00004009 ; 1513.8641730465356
   & &BD3BA74E
   & &3B9E90C8
   & &00004006 ; 233.17823205143102
   & &E92DA09D
   & &9FFECB3E
   & &0000400B ; 4368.0886700674173
   & &8880B598
   & &A6FD8544

       
exm1 STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 BIC R2,R0,#&80000000
 SUB R2,R2,#&3F000000
 LDFE F1,LOG2E
 CMP R2,#&00D00000
 BGE EXM1L
 MUFE F0,F0,F1
exm1y
 ADR R3,EXM1COFFS
 MUFE F1,F0,F0
; CALCULATE P
 lfm F2,2,[R3],#24   ; P2,P1
 MUFE F2,F2,F1
 ADFE F2,F2,F3
 MUFE F2,F2,F1
 lfm F3,1,[R3],#12   ; P0
 ADFE F2,F2,F3
 MUFE F0,F0,F2       ; x * P
 lfm F2,2,[R3],#24    ;Q1 Q0
 ADFE F2,F1,F2
 MUFE F2,F2,F1
 ADFE F2,F2,F3    ;  Q
 SUFE F1,F2,F0    ; Q - xP
 MUFE F0,F0,#2    ; 2 x P
 DVFE F0,F0,F1
 mov pc, r14

EXM1L STMFD SP!,{R14}
 BL yexpf
 SUFE F0,F0,#1
 LDMFD SP!,{PC}

;      ----------

ypow TEQ R0,#0
 BMI POWNEG
 TEQEQ R1,#0
 BEQ POWZER
 STMFD SP!,{R2,R3,R14}
 LDFD F0,[SP]
 URDE F1,F0
 NRME F1,F1
 SUFE F1,F0,F1
 CMF F1,#0
 BNE POWUSU
 ABSE F1,F0
 LDFS F2,POTRESH
 CMF F1,F2
 BGT POWUSU

 FIX R2,F0
 STMIA SP,{R0,R1}
 LDFD F1,[SP],#12   ; WE KILL R14 ON STACK AS WELL
 TEQ R2,#0
 BEQ PPZ
 MVFS F0,#1
 TEQ R2,#0
 MOVEQS PC,R14
 RDFMIE F1,F1,#1
 RSBMI R2,R2,#0

POWLOOP MOVS R2,R2,LSR #1
 MUFCSE F0,F0,F1
 MOVEQS PC,R14
 MUFE F1,F1,F1
 B POWLOOP

PPZ CMF F1,#0
 ADFVSE F0,F1,#0  ; RETURN IF nan
 MOVVSS PC,R14
 SUB IP,R0,#&0FF00000
 TEQ IP,#&70000000
 BEQ PAWX           ; +- INF ** 0
 MVFSZ F0,#1
 mov pc, r14


POWNEG
 STMFD SP!,{R0,R1}
 LDFD F2,[SP],#8
 CMF F2,#0
 BVS PNAN
 MOVS IP,R2
 RDFMIE F0,F2,#1
 MVFPLE F0,F2
 BICMIS IP,IP,#&80000000
 TEQEQ R3,#0
 BEQ PMZ
 MVFez F1,F0
 MOV R0,IP,LSL #12
 MOV IP,IP,LSR #20
 SUB IP,IP,#&FF
 TEQ IP,#&700
 BEQ PF1            ;  INF OR NAN
 SUBS IP,IP,#&300
 BMI PF1            ; .. fraction
 ORR R0,R0,R3,LSR #20
 MOV R1,R3,LSL #12
 BEQ PN2
PN1 MUFE F0,F0,F0
 TEQ R0,#0
 MUFMIE F0,F0,F1
 MOVS R1,R1,LSL #1
 ADCS R0,R0,R0
 SUBS IP,IP,#1
 BNE PN1
PN2 TEQ R0,#0
 TEQEQ R1,#0
 MOVEQS PC,R14

PF1 MVFez F0,F2      ; original first arg
 STMFD SP!,{R2,R3}
 LDFD F1,[SP],#8 ; original scnd arg
 POWD F0,F0,F1
 mov pc, r14

PMZ ;  handle (-x)^0
 BICS IP,R0,#&80000000
 TEQEQ R1,#0
 BEQ PAWX      ; (-0)^0
 MOVS IP,IP,LSR #20
 SUBS IP,IP,#&FF
 TEQ IP,#&700
 BEQ PAWX
 MVFsz F0,#1
 mov pc, r14

POWZER ;  handle 0 ^ X
 STMFD SP!,{R2,R3}
 LDFD F0,[SP],#8
 CMF F0,#0
 ADFVSE F0,F0,#1
 MOVVS PC,R14
 BEQ PAWX
 MVFsz F0,#0
 RDFMIE F0,F0,#1
 mov pc, r14

PAWX ; supply an invalid operation
 MNFsz F0,#1
 LOGD F0,F0
 mov pc, r14

PNAN STMFD SP!,{R2,R3}
 LDFD F1,[SP],#8
 ADFE F0,F1,F2
 mov pc, r14

POWUSU BL ylog
 LDFD F1,[SP],#8
 MUFE F0,F0,F1
 LDMFD SP!,{R14}
 B yexpf

POTRESH DCFS 65000.0




twotox STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
twotoxf LDFS F2,TOVT
 ABSEZ F1,F0
 CMF F1,F2
 ADR R3,P2
 BMI YEXPZ
 B OVER_UNDER


yexp
 STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
yexpf LDFS F2,OVT
 ABSE F3,F0
 ADR R3,LOG2E
 lfm F1,1,[R3],#12
 CMF F3,F2
 MUFE F0,F0,F1
 BPL OVER_UNDER
YEXPZ URDE F1,F0
 NRME F1,F1    

 FIX R2,F1

 SUFE F0,F0,F1
; -1/2 <= F0 <= 1/2
 MUFE F1,F0,F0
; CALCULATE P
 lfm F2,2,[R3],#24   ; P2,P1
 MUFE F2,F2,F1
 ADFE F2,F2,F3
 MUFE F2,F2,F1
 lfm F3,1,[R3],#12   ; P0
 ADFE F2,F2,F3
 MUFE F0,F0,F2
 lfm F2,2,[R3],#24    ;Q1 Q0
 ADFE F2,F1,F2
 MUFE F2,F2,F1
 ADFE F2,F2,F3
 SUFE F1,F2,F0
 ADFE F0,F0,F2
 DVFE F0,F0,F1
 TEQ R2,#0


 MOVEQS PC,R14
 CMP R2,#&4000
 BGE OVERFL
 ADD R0,R2,#&FF
 ADDS R0,R0,#&3F00
 MOV R1,#&80000000
 MOV R2,#0
 BLE DENORM_OR_ZER
EX2 STMFD SP!,{R0,R1,R2}
 lfm F3,1,[SP],#12
 MUFE F0,F0,F3
 mov pc, r14
 
TOVT DCFS 16448.0
OVT DCFS 11400.0
LOG2E & &00003FFF
      & &B8AA3B29
      & &5C17F0BC
P2    & &00003FF9
      & &BD2E42AB
      & &70BDAA7B
P1    & &00004003
      & &A19DD498
      & &9F60DB06
P0    & &00004009
      & &BD3D047F
      & &734DBD67
Q1    & &00004006
      & &E92F287A
      & &E89542C1
Q0    & &0000400B
      & &8881B17C
      & &3A652AD2


DENORM_OR_ZER
 RSB IP,R0,#1
 MOV R3,#&80000000
 MOV R1,R3,LSR IP
 SUBS IP,IP,#32
 MOVPL R2,R3,LSR IP
 B EX2

OVER_UNDER
 BVS ENAN
 CMF F0,#0
 MVFLTE F0,#0
 MOVLTS PC,R14
 ; we would need 2**f0 at this point but we have an overflow so we dont care
 EXPD F0,F0
 mov pc, r14

OVERFL
 LDFS F0,TOVT
 EXPD F0,F0
 mov pc, r14

ENAN ADFE F0,F0,#0
 mov pc, r14

; -----------------------------------------------------------


ylog10 STMFD SP!,{R14}
       BL ylog
       LDFE F1,L10E
       MUFE F0,F0,F1
       LDMFD SP!,{PC}
L10E & &00003FFD
     & &DE5BD8A9
     & &37287195

ylog
 LDFD F1,SQX2
 MOVS R2,R0,ASR #20
 BMI LSUP
 MOV R3,#&FF
 ORR R3,R3,#&300
 BEQ LGZER_OR_DEN
 BIC R0,R0,R2,LSL #20
 ORR R0,R0,R3,LSL #20
 SUB R2,R2,R3           
 TEQ R2,#&400          ; INFINITY OR NAN AS INPUT
 BEQ LSUP1
LN1 STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8

; F0 IS NOW BETWEEN 1 AND 2
; F1 CONTAINS SQRT2
; R2 IS ORIGINAL EXPONENT
 CMF F0,F1
 MUFGTE F0,F0,#0.5
 ADDGT R2,R2,#1
; WE NEED NOW F0 = (F0-1)/(1+F0)
 SUFE F1,F0,#1
 ADFE F2,F0,#1
 ADR R1,LGCO7
 DVFE F0,F1,F2
 lfm F2,2,[R1],#24
 MUFE F1,F0,F0  
 MOV R0,#5
 MUFE F2,F2,F1
 ADFE F2,F2,F3
LG2
 MUFE F2,F2,F1
 lfm F3,1,[R1],#12
 SUBS R0,R0,#1
 ADFE F2,F2,F3
 BPL LG2
 TEQ R2,#0
 MUFE F0,F0,F2
 FLTNEE F1,R2
 MOVEQS PC,R14
 LFM F3,1,[R1],#12
 MUFE F1,F1,F3
 ADFE F0,F0,F1
 mov pc, r14

SQX2 DCFD 1.41421356237309504880
LGCO7    & &00003FFC
       & &98AA3FC1
       & &56F198AE
ZCO6    & &00003FFC
       & &9CB03E62
       & &9FC48DB6
       & &00003FFC
       & &BA34901A
       & &A65F3672
       & &00003FFC
       & &E38E20BB
       & &67572CF5
       & &00003FFD
       & &924924AD
       & &BBFE1303
       & &00003FFD
       & &CCCCCCCC
       & &ABC1FB26
       & &00003FFE
       & &AAAAAAAA
       & &AAB39730
       & &00003FFF
       & &FFFFFFFF
       & &FFFFFF4C
LOG2   & &00003FFE
       & &B17217F7
       & &D1CF79AC


LGZER_OR_DEN ;
 TEQ R0,#0
 TEQEQ R1,#0
 MVFEQsz F0,#0
 BEQ LSP1      ; input was zero
 RSB R2,R3,#1
LGNORM MOVS R1,R1,LSL #1
 ADC R0,R0,R0
 SUB R2,R2,#1
 TST R0,#&100000
 BEQ LGNORM
LGNORM1
 ORR R0,R0,R3,LSL #20
 B LN1

LSUP1 ADD R0,R0,#&40000000
LSUP STMFD SP!,{R0,R1}
     LDFD F0,[SP],#8
LSP1 LGNE F0,F0
     mov pc, r14


ycos STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 ADR R1,sCO8
 lfm F1,1,[R1,#-12]
 MUFE F0,F1,F0
 RSFE F0,F0,#0.5
 B YSIN1

ysin STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 ADR R1,sCO8
 lfm F1,1,[R1,#-12]
 MUFE F0,F1,F0
YSIN1 MOV R0,R0,LSL #1
 CMP R0,#&83000000
 URDE F1,F0
 MOV R0,#7
 NRME F1,F1
 BHI INVARG
 FIX R2,F1
 SUFE F0,F0,F1
SIN1 TST R2,#1
 MUFE F1,F0,F0
 lfm F2,1,[R1],#12
 MNFNEez F0,F0
SIN2
 MUFE F2,F2,F1
 lfm F3,1,[R1],#12
 SUBS R0,R0,#1
 ADFE F2,F2,F3
 BPL SIN2
 MUFE F0,F0,F2
 mov pc, r14

INVARG 
 CMF F0,#0
 BVS SINAN
 SUFE F0,F0,F1   
 MUFE F1,F1,#0.5    ; = -1  -.5  0  .5  1
 URDE F2,F1
 NRME F2,F2         ; = -1  -1/0 0  0/1  1
 SUFE F1,F1,F2      ; = 0   +-.5 0  +-.5 0
 MUFE F1,F1,#2      ; = 0   +-1  0  +-1  0
 FIX R2,F1
 B SIN1


OPI    & &00003FFD   ; 1/PI
       & &A2F9836E
       & &4E44152A
sCO8 & &00003FEA
    & &CE789708
    & &343B50C6
sCO7 & &80003FEF
    & &B7BC195F
    & &5A69B150
sCO6 & &00003FF3
    & &F479AB60
    & &12957FD1
sCO5 & &80003FF7
    & &F183A6E1
    & &E13B382E
sCO4 & &00003FFB
    & &A83C1A42
    & &75F16822
sCO3 & &80003FFE
    & &99696673
    & &137E85C4
sCO2 & &00004000
    & &A335E33B
    & &AD531C80
sCO1 & &80004001
    & &A55DE731
    & &2DF290FD
sCO0 & &00004000
    & &C90FDAA2
    & &2168C231

SINAN ADFE F0,F0,#0
 mov pc, r14
; ------------------
ytan STMFD SP!,{R0,R1}
 LDFD F0,[SP],#8
 LDFE F1,OPI
 MUFE F0,F1,F0
 URDE F1,F0
 NRME F1,F1
 SUFE F0,F0,F1
 LDFS F1,ONE_O_FOUR
 MOV R2,#0
 CMF F0,F1
 RSFGTE F0,F0,#0.5
 MOVGT R2,#1
 CNF F0,F1
 MNFLTE F0,F0
 SUFLTE F0,F0,#0.5 
 MOVLT R2,#1
 MUFE F1,F0,F0
 ADR R1,TANCOFFS   ; 4286  P3 P2 P1 P0 Q3 Q2 Q1 Q0
 LDFE F2,[R1],#12
 MOV R0,#2
TAN2
 MUFE F2,F2,F1
 LDFE F3,[R1],#12
 SUBS R0,R0,#1
 ADFE F2,F2,F3
 BPL TAN2
 MUFE F0,F0,F2   ; x*P
 LDFE F2,[R1],#12
 MOV R0,#2
 ADFE F2,F1,F2   ; x^2 + q(3)
TAN3
 MUFE F2,F2,F1
 LDFE F3,[R1],#12
 SUBS R0,R0,#1
 ADFE F2,F2,F3
 BPL TAN3
 TST R2,#1
 DVFEQE F0,F0,F2
 DVFNEE F0,F2,F0
 mov pc, r14


ONE_O_FOUR DCFS 0.25
TANCOFFS
    & &80004002 ; -2   -45.649319438665628
    & &B698E731
    & &E6E763D6
    & &00004006 ; -6   14189.854252761779
    & &DDB76AC1
    & &3C696C66
    & &80004008 ; -a -895888.44006768044
    & &DAB9070A
    & &84687642
    & &00004008 ; -e   10888600.437281687
    & &A625986F
    & &F1B1530C
    & &80004004 ; -4  -1014.6561902528854
    & &FDA9FF05
    & &67062D0D
    & &00004008 ; -8    135382.71280511908
    & &8435AD9E
    & &995CBB7B
    & &80004008 ; -c   -3991309.518035165
    & &F39C3612
    & &77CF7200
    & &00004006 ; -10    13863796.663567629
    & &D38B74A9
    & &DF9171F3


 END
